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Abstract 

Precise variational solutions are given for problems involving diverse fermionic 
and bosonic = 2— 7-body systems. The trial wave functions are chosen to be 
combinations of correlated Gaussians, which are constructed from products of 
the single-particle Gaussian wave packets through an integral transformation, 
thereby facilitating fully analytical calculations of the matrix elements. The 
nonlinear parameters of the trial function are chosen by a stochastic technique. 
The method has proved very efficient, virtually exact, and it seems feasible for 
any few-body bound-state problems emerging in nuclear or atomic physics. 
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I. INTRODUCTION 



Few-body problems of interacting particles have vital importance in all branches of 
physics from hadron to celestial levels. The main interest in the few-body problems lies 
in, e.g., finding an accurate solution for the system, testing the equation of motion and 
the conservation laws and symmetries, or looking for unknown interactions governing the 
system. 

The investigation of few-nucleon systems interacting via realistic forces has always been 
in the center of the interest. Considerable effort has been exerted to obtain accurate ground- 
state properties of the few-nucleon systems with Faddeev-Yakubovsky (FY) variational 
variational Monte Carlo (VMC) 0J| and Green's function Monte Carlo (GFMC) 
methods 0,0. Most of these approaches has focused on three- or four-body problems. 

To treat an A^-particle system, one needs to cope with a large number of variables required 
to specify the wave function. By using (N—l) relative coordinates to describe the system, for 
example, the discretization on a mesh with p points, or the expansion of the function of the 
relative motion between the particles in terms of p suitably chosen functions leads to 
mesh points or basis functions, which becomes prohibitively large with increasing N. All 
but the Monte Carlo methods face this difficulty as the number of particles increases. The 
VMC and GFMC methods have proved to be most successful by being able to go beyond 
the four-nucleon problem The secret of the efficiency of the Monte Carlo methods is 



the use of an importance sampling of the most relevant parts of the configuration space. 
This fact naturally raises a question: Even if the wave function of the iV-particle system is 
expanded into an (excessively) large number of basis functions, can't one reduce the problem 
to a tractable one by selecting "the most important" basis functions ? 

The aim of this paper is to present an alternative variational approach, the stochas- 
tic variational method (SVM) [|12|,|13[, by using the correlated Gaussians as basis functions 



1J,|15[. Examples whose solutions were known before are used to demonstrate the perfor- 
mance of the method in treating nuclear as well as Coulomb interactions. To highlight some 
new physics, we have also included problems that have been hitherto unsolved. We give 
the formulation and some details of the method of the calculation and show applications to 
N = 2 — 7-particle systems. 

The variational foundation for the time-independent Schrodinger equation provides a 
solid and arbitrarily improvable framework for the solution of bound-state problems. The 
crucial point of the variational approach is the choice of the trial function. There are two 
widely applied strategies: (1) to select the most appropriate functional form to describe the 
short-range as well as long-range correlations and to compute the matrix elements by Monte 
Carlo technique, or (2) to use a number, possibly a great number, of simple terms, which 
facilitate the analytical calculation of the matrix elements. We follow the second course by 
using an expansion over a correlated Gaussian "basis" . 

To solve the A^-particle problem, it is of prime importance to describe the correlation 
between the particles properly. The correlation is conveniently represented by a correla- 



tion factor, F = Yli^j fij [^|5|,|7|,|8|,p!6[] . Most calculations have used this form of F directly 
to evaluate the matrix elements. Such calculations are, however, fairly involved beyond 
the three-particle system and performed by Monte Carlo integrations. An alternative way 
to incorporate the correlation is to approximate fij as a linear combination of Gaussians 
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exp{—aij{ri — rj)"^). The A^-particle basis function then contains product of these Gaus- 
sians: ni<j exp(— ajj(rj — rj)"^) = exp(— X^i^j ctijlrj — ^j)"^)- These Gaussian functions are 



widely used in variational calculations (see, for example, We will apply a more 

general form of the correlated Gaussian functions which allow for nonzero orbital angular 
momentum, and will use the more convenient Jacobi relative coordinates instead of the 
relative distance vectors. The correlated Gaussians have an important advantage. Their 
Hamiltonian matrix elements can be analytically calculated in a unified framework, thus 
enabling one to avoid the formidable calculation involving the correlation factor F. 

The variational approximation, however, may run into difficulties for the following rea- 
sons: (i) if the nonlinear parameters specifying the basis functions are varied, it is difficult 
to optimize them, (ii) if they are not, then the number of terms required may be excessively 
large, and, in both cases, (iii) the trial function of proper symmetry becomes extremely 
involved. For example, conventional methods |]6|JT5[1 for the choice of the Gaussian parame- 
ters lead to prohibitively large bases for more than 3 or 4 particles, which has limited the 
applicability of the Gaussian basis to few-body problems. 

One can circumvent the optimization problem including large number of nonlinear pa- 
rameters or the diagonalization of huge matrices by using the SVM. The SVM attempts 
to set up the most appropriate basis functions by the following stepwise procedure: One 
generates a would-be basis function by choosing the nonlinear parameters randomly, judges 
its utility by the energy gained by including it in the basis, and either keeps or discards it. 
One repeats this "trial and error" procedure until the basis set up leads to convergence. The 



original procedure of the SVM, proposed in [|1^], has recently been developed further and 



successfully applied to multicluster descriptions of light exotic nuclei, such as 'Re=a + n + n, 
®He=a + n + n + n + n, ^Li=a + t + n + n, and ^C=a+^B.e+p + p ||T3|JT^ . Learning from 
these applications, we have now generalized and refined the method further to encompass 
diverse systems emerging in nuclear and atomic physics. 

Besides the large number of nonlinear parameters, the treatment of the increasing num- 
ber of partial waves in the expansion of the wave function would also pose a formidable 
task. We propose here an alternative formulation to cope with this problem. Instead of us- 
ing the partial wave expansion, the angular dependence of the wave function is represented 
by a single solid spherical harmonics whose argument contains additional variational param- 
eters. This form makes the calculation of the matrix elements for nonzero orbital angular 
momentum much simpler than other methods. 

It will be demonstrated that the present method has several unique features: It is based 
on a fully analytical calculation for most types of interactions and thus ensures high accuracy 
and speed. Its calculational scheme is quite universal and needs no change depending on 
whether the system contains nuclear or Coulombic or other interactions. It has no difficulty 
in treating the system of particles of unequal masses. More importantly, the wave function 
is obtained in a compact, analytical form and thereby can be readily used in calculations of 
physical properties. 

As you will see later, the present method has turned out to be very accurate, and we think 
it is worth while to make the method and the results easily available and reproducible for 
interested readers. We collect all the needed ingredients of our method in order. Some of the 
formulae are our original developments or generalizations of known relations to iV-particle 
matrix elements, and some others are collected here to make the paper self-contained. The 
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calculation of the matrix elements presented here is different from the one of Refs. ||I8|,[T9 



m 



many aspects: The motion of the centre-of-mass is removed from both the Hamiltonian and 
the wave function. Two-particle potential matrix elements of arbitrary radial form factor 
are evaluated in a unified way by reducing them to the calculation of appropriate correla- 
tion functions corresponding to the interaction. The calculation of the matrix elements is 
extended to nonzero orbital angular momentum as well. The symmetrization postulate is 
imposed on the wave function at the single-particle level, which provides several advantages, 
especially in evaluating the matrix elements of state-dependent realistic nuclear interactions. 

The organization of the paper is as follows. Section II defines the correlated Gaussian 
basis functions and gives the details of the stochastic procedure of selecting the basis set. 
Section III contains the method of calculating the matrix elements. The main steps are the 
calculation of matrix elements in Slater determinants (or permanents for bosons) consisting 
of single-particle Gaussian wave packets, the elimination of the center-of-mass motion with 
a very simple manipulation, and the transformation to the correlated Gaussian basis. This 
section also presents the modifications needed for treating systems of particles of unequal 
masses. Section IV presents numerical results for various systems of particles which interact 
via nuclear potentials or power-law potentials. Section V gives a brief summary. In the 
appendices the most important auxiliary formulae are collected to facilitate any future use 
of the formulation. 



II. THE CORRELATED GAUSSIANS AND THE STOCHASTIC VARIATIONAL 

METHOD 

A. Basis functions 

Since the variational method is always limited by the form chosen as a trial function, the 
trial function must be flexible enough to be able to describe the full variety of correlations 
between the nucleons, e.g., the short-range correlation due to the strong repulsive force, the 
a-clustering typical in some light nuclei, or the long-range correlation at large distances in 
light halo nuclei. The correlation between the nucleons can be described by functions of 
appropriate relative coordinates. 

Any square-integrable function with angular momentum Im can be approximated, to any 
desired accuracy, by a linear combination of nodeless harmonic-oscillator functions (Gaus- 
sians) of continuous size parameter a: 

r,^(r) ~e-^'^'^'3^i^(r), with 3^;^(r) = r'y^^(r). (1) 

A generalization of this to A^-nucleon systems contains a product of the Gaussians as men- 
tioned in the previous section. It is convenient to use a set of the Jacobi coordinates 
X = (xi, ...,X7v_i), instead of N{N—l)/2 relative distance vectors (r^ — r^). An A^-nucleon 
basis function, a so-called correlated Gaussian, then looks like 

'^{ls)jmtmt{^, A) = ^{e-^*^^ [0l{^)Xs]jm ^tmt}. (2) 

where x, the transpose of x, stands for the row vector comprising the Jacobi coordinates, x 
and X are the spin and isospin functions. A is an [N—l)x[N—l) positive-definite, symmetric 
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matrix of nonlinear parameters, specific to each basis element, and the quadratic form, ^Ayi, 
involves scalar products of the Cartesian vectors: 

N-l N-1 

i^x = ^ ^ AijXi ■ Xj. (3) 
i=i j=i 

The operator A is an antisymmetrizer defined by 

'4=^E«ig^(^)^' (4) 

where the sum runs over all permutations of the N nucleon indices and sign(P) stands for 
the parity of the permutation P. For a system of identical bosons, the antisymmetrizer is 
to be replaced with a symmetrizer. For a general case A is to represent the operator that 
imposes the proper symmetry on the wave function. 

The function 6lml{'^) in Eq- (2), which represents the angular part of the wave function, 
is a generalization of y and can be chosen as a vector-coupled product of solid spherical 
harmonics of the Jacobi coordinates 

Each relative motion has a definite angular momentum in Eq. (5). It may be important, 
however, to include several sets of angular momenta {li, I2, ■■■■ In-i] -^12, -^'1235 •••) for a realis- 
tic description. The various possible partial wave contributions increase the basis dimension; 
moreover, the calculation of matrix elements for this choice of 6lml{'^) becomes too compli- 
cated. This choice is apparently inconvenient especially as the number of nucleons increases. 
To avoid this, we propose a different choice as the generalization of y-. 

N-l 

OlmA^) ^Vklml{u,x) ^v^^'^^YlmAv), with v = ^ li^Xj. (6) 

1=1 

Only the total orbital angular momentum appears in this expression and it contains a 
parameter u = {ui, ...,un^i). The vector u may be considered as a variational parameter 
and one may try to minimize the energy functional with respect to it. It defines a linear 
combination of the Jacobi coordinates, v, and the wave function of the system is expanded 
in terms of its angle v. The minimization amounts to finding the most suitable angle or a 
linear combination of angles. The factor of v'^^'^^ plays an important role in improving the 
short-range behavior of the wave function. A remarkable advantage of this form of 9lml{^ 
is that the calculation of matrix elements becomes much simpler than in the former case 
because the coupling of (A^— 1) angular momenta is completely avoided. 

The two forms of 9lml{^) f^-ct closely related to each other. Any of the functions 

of Eq. (5) may be expressed in terms of a linear combination of the terms, v'^^'^^Y^Ml ( v) , 
by using some appropriate sets of u values provided that each term satisfies the condition 
2K + L<li+ ■ ■ ■ +In-i and contains a monomial of degree h + ■ ■ ■ + In-i — '^K — L in the 
variables, xi^, ...,XAr_i^. Therefore, if one can calculate the matrix elements using 9lml{^) 
defined in Eq. (6), then those with the previous form of 9lml{'^) can be obtained readily. 
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The correlated Gaussian basis with the function 6'la/^(x) of Eq. (6) has parity (— 1)'''. 
To construct a function with parity (— l)^"*"^, Eq. (6) must be shghtly generahzed, e.g., to 

6'LAfz,(x) = [r7i^L(u,x)r7oi(u',x)]iAfz,- (7) 

To assure positive definiteness, the matrix A in Eq. (2) is in general expressed as A = 
GA'G, where G is an (A^— l)x(A^— 1) orthogonal matrix containing (A^— l)(iV— 2)/2 parameters 
and A' is a diagonal matrix, {A')ij = a'-6ij, including (A^— 1) positive parameters a'-. Although 
no restriction on the parameters of the matrix G is in principle necessary, it is advisable to 
avoid too many variables if possible. The most naive choice would be to take G as a unit 
matrix, which is equivalent to using only a single set of the Jacobi coordinates, and then 
to try to reach convergence by including higher partial waves successively. Many examples 
show 1^,17 1, however, that this does not work well because the convergence is generally slow 



and moreover the computational cost of using high partial waves is quite expensive. 

The matrix G can also be chosen as one of the rotation matrices that connect the set of 
the Jacobi coordinates to other sets of independent relative coordinates. Figures la-lf show 
all topologically different sets of independent relative coordinates for a system of six identical 
particles. The set of coordinates in Fig. la is what we call the set of the Jacobi coordinates 
x. A correlation conforming to a specific set of relative coordinates x' = (x'^^, ...,x'^_^) can 
be most efficiently described by tailoring the form of the basis function to this set of relative 
coordinates, that is, by using the form, exp{ — | Ylj^^ oli^i "x'j}. Since the coordinates x' can 
be obtained by an appropriate rotation IZ of the Jacobi coordinates as x' = T^x, the basis 
function of such type can be clearly encompassed in the trial function of Eq. (2) by choosing 
G = TZ. The correlated Gaussian basis thereby allows for various correlations between the 
nucleons and different asymptotics at large distances flexibly. Depending on the character 
of the problem a more general choice of G might be necessary. 

By selecting a set of basis functions [ipi] i = 1, /C} [^Ji = '4'(LiSi)JMTMTi^y ^i)] that ad- 
equately spans the state space, the wave function of the A^-nucleon system can be expanded 
as 

ic 

^ = Ec.^^' (8) 

i=l 

where c = (ci, c^c) is the set of linear variational parameters. The Ritz variational method 
defined by this trial function reduces to the generalized algebraic eigenvalue problem 

He = EAfc, (9) 

where Ti and Af are, respectively, the matrices of the Hamiltonian and of the overlap 

nij = {^i\H\i;,), and Mj = (Hi^j) {t, J = I, IC). (10) 



B. Stochastic selection of parameters and solution of the eigenvalue problem 

Since a linear combination of the correlated Gaussians forms a dense set, there are 
different sets of A that represent the wave function equally well. This enables one to select 
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the most appropriate parameters randomly. We set up the basis stepwise by choosing A 
from a preset domain of the parameter space and increase the basis dimension one by one. 
In the first step we select a number of parameter sets A randomly, and keep the one that 
gives the lowest energy. Next we generate a new random set and calculate the energy with 
this two-element basis. As one more basis state always lowers the energy, we quantify its 
"utility" by the energy gained by including it in the basis. If the energy gain is larger 
than a preset value, e, then we admit this state to the basis, otherwise we discard it and 
try a new random candidate. This is repeated until the energy converges. The rate of 
convergence can be controlled by dynamically decreasing the value of e during the search. 



This procedure is more advantageous than the earlier versions []T2|,[T3| and, although not a 
full optimization, results in very good and relatively small bases. A similar procedure, called 
"stochastic diagonalization" has been used to determine the smallest eigenvalue of extremely 
large matrices . 

To have an economical algorithm for setting up the basis by a trial and error method, one 
has to find an efficient way to solve the eigenvalue problem, Eq. (9). The full diagonalization 
is rather time consuming and in fact unnecessary because (i) in the (/C+l)th step of the 
procedure we can use the result of the /Cth step and (ii) to judge the usefulness of a would- 
be basis state, only the lowest eigenvalue is needed. Let us assume that in the /Cth step 
the Hamiltonian matrix Ti. is diagonalized; its eigenenergies are Ei < E2 < ... < E/c and 
its corresponding normalized eigenfunctions are The eigenvalue problem in the 

(/C-|-l)th step takes the form 

/ Et ... 



Ex 




\{llj!C+i\H\^l) 



E, 



X 



{^!c\H\^,c+i) 
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... {i)^+x\ii\^K) {ll)K+l\H\i)K+l) ) \cK+ll 
... (^i|^yc+i) \ / ci \ 







^1) 



1 (^^l^^+i) 

... {i)K+l\^K) (^/C+l|^/C+l) ) 



V C/C+1 J 



(11) 



By using the Gram-Schmidt orthogonalization method, that is, by defining 



IV'^+i)-Ef=i 






Vl/,)(v[/,|V.^+l))l/2 



(12) 



this generalized eigenvalue equation can be reduced to the conventional form 



/ El ... gi \ / ci \ 





V ^1 



••• Ejc qjc 
... Qk a ) 



V Cjc+i J 



E 



(13) 
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where 



= {<i/i\H\^,c+i), a = {i>,c+i\H\^ic+i)- (14) 
The eigenvalues are easily obtained by finding the roots of the secular equation 

\{E) ^ n(i?^ -E)Ua-E)-Y.^ -^^\ = 0. (15) 

This secular equation has (/C+1) roots {E'fji = 1,...,/C+1} fulfilling the inequalities E'l < 
El < E'2 < E2 < ■■■ < Etc < E'lc+i- The eigenvectors are readily obtained after substituting 
the eigenvalues {E'i]i = 1, ...,/C+l} into Eq. (13). Note that one has to determine only the 
lowest eigenvalue E'l for the admittance criterion. 

III. CALCULATION OF THE MATRIX ELEMENTS 

In this section we will give the details of the method of calculating the matrix elements 
between the basis function of Eq. (2). The calculation consists of three steps: (A)The 
calculation of the matrix elements between the Slater determinants of the Gaussian wave- 
packet single-particle functions, (B)A transformation from the single-particle coordinate 
representation to the relative and center-of-mass coordinate representation, (C)An integral 
transformation from the Gaussian wave-packet functions to the correlated Gaussian basis. 
A procedure similar to steps (A) and (B) was used to manipulate algebraically the antisym- 
metrization operation and the transformation of the coordinates for complex cluster systems 



21 1 . In step (A) the Slater determinant for the A^-nucleon wave function is constructed by 
distributing the nucleons at positions (si, ...,S7v). These position vectors serve as the gen- 
erator coordinates. The Slater determinant of the Gaussian wave packets is often used in 



nuclear theory, e.g., in cluster model pT|-p^ and fermionic or antisymmetrized molecular 



dynamics [^,0 . The Hamiltonian matrix elements are analytically evaluated with the use 
of technique of the Slater determinants P7| , p2| , and can be expressed as a function of the 



generator coordinates. In step (B) the center-of-mass motion is completely separated from 
the intrinsic motion, and thus the trial wave function acquires the translational invariance. 
The separation of the center-of-mass motion is particularly simple in this formulation. In the 
last step (C) the matrix elements expressed in terms of the intrinsic generator coordinates 
are transformed to those between the correlated Gaussian basis functions with a definite 
angular momentum. Some of the essential parts of the calculational scheme is our original 
development, and some of them is a generalization of the technique used in the nuclear 
cluster model (see, for example, [^). We also show in subsection III.D those modifications 



which are needed to treat the system of particles of unequal masses. 

A. Slater determinants of Gaussian wave packets 

The ith nucleon with mass m, spin cTj and isospin is to be put in the single-particle 
Gaussian wave packet 
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Vl.^niji) = <(ri)xi.,'^i.,, (16) 

with 

<(r.)-{^fe-^i--r, aud . = (17) 

where is the position vector of the nucleon, xi^. and are its spin and isospin function. 
The angular frequency u is not a variational parameter and may be taken an arbitrary con- 
stant. The Si parameter or "generator" coordinate will be used in an integral transformation 
to derive the matrix elements between the Gaussian basis functions. A Slater determinant 
of these Gaussian packets is defined by 

</)«(si,...,s^) = ^|n^:^...,(r.)|, (18) 

where n = (ctiTi, ...,a]\fT]\r) is the set of the spin-isospin quantum numbers of the nucleons. 
The spins and the isospins of the nucleons are successively coupled to add up, respectively, 
to the total spin SMs and isospin TMj- of the A^-nucleon system: 



XSMs 



SMs 



TMt 



To simplify the notation, the intermediate quantum numbers are suppressed in the following. 
The wave function in the "generator coordinate space" with the definite spin and isospin 
quantum numbers is a linear combination of the Slater determinants of the Gaussian packets: 

^SMsTT, (Sl, Sat) = A{(p''^^{ri)...(p''^^{rN)XSMs^TMT} = ^^^/.(Sl, Sat), (20) 

K 

where is a product of the Clebsch-Gordan coefficients needed to couple the spin and 
isospin as defined in Eq. (19). 

The Hamiltonian of the A^-nucleon system reads as 

N 2 N 

^ = E|j- + E% (21) 

i=i i<j 

The matrix elements of the Slater determinants can easily be evaluated using the well-known 
rules [0,0 . To make this paper self-contained, we have collected all the needed ingredients 
in Appendices A, B and C. The overlap of the Slater determinants is found to take the form 



1 (o) 

{^SMsTmA^U ...,Sn)\^SMsTMt{^'i, -,8^)) = ^cf^c'^^^ ^ (22) 



G'>-e-2=^ 

i=l 

where A^f* is a 2A^ x 2A^ real, symmetric matrix and s stands for the 2A^-dimensional row 
vector comprising the single-particle generator coordinates, (si, s^r, s'^, s'^). To simplify 
the notation, we refer to the set of the vectors (s'^^, s'^) alternatively as (sjv+i, S2Ar). 
Note, therefore, that the quadratic form, sA-°''s, reads as 
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2N 2N 

M°^s = EEKO,,«.-«^- (23) 

j=i k=i 

The matrix elements of the kinetic energy operator can also be expressed in terms of the 
same 's and A\"^ 's as 

^ p2 

{^SMsTMri^l^ Sat)! E | ^5MsTMt(s'i, s'iv)) 

i=l 

= ^tci''>{lN-\,Ar>s)e-i'-r>.. (24) 

The matrix elements of any term of the two-body interaction can be expressed as an integral 
of the two-particle correlation function multiplied by the radial form factor, V{r), of the term 
Vij as below: 



N 

($SMsTMt(Si, ■■■,SN)\Y,^ij\^S'M'gT'M^{s'i, s'jv)) 



rip 

= / drV{r)e-'''" C^^^Piis, r)e-5^^i''^+'i'-^ (25) 

•' i=l 

where Pi{s, r) is a polynomial of s and r, A^^ a 2Nx2N symmetric matrix, and dj • r takes 
the form 

27V 

di-r = ^P(i)jSj-r. (26) 

3=1 

The polynomial part reduces to unity (Pj(s, r) = 1) in the case of pure central forces, 
but it has a rather simple form for spin-orbit, tensor, and other interactions as well. See 
Appendices B and C. The Cf''s, A^^'s, Cf^^'s, and ^^^^'s, etc. are obtained with the use of 
mathematical manipulation languages or fortran programs. See Appendix C. 



B. Transformation to relative and center-of-mass coordinates 

To eliminate the center-of-mass motion, we transform the single-particle coordinates to 
the relative and center-of-mass coordinates. For this purpose we choose one particular set of 
relative coordinates, the Jacobi coordinates, which is expressed in terms of the single-particle 
coordinates as 

N 

^i = Y.UikVk (^ = 1,...,A^), (27) 

fc=i 

where the transformation matrix U is defined by 
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u 



... \ 



■i 1 ... 



1 



V 



AT-l N-l 
1_ 1_ 

N N 



N 



I 2 



and U 



1 

2 

^ 

" 3 



3 - 1\ 

i _i 1 

3 Af 



V 



(28) 



1/ 



Similarly, the single-particle generator coordinates are transformed to the relative and center- 
of-mass generator coordinates: 



AT 



N 



(29) 



fe=i 



k=l 



The reduced masses corresponding to the transformation U are given by 



i + 1 



-m 



l,...,N—l), and /in = Nm. 



(30) 



The product of the Gaussian single-particle wave packets can then be written as a product 
of Gaussians depending on the relative and center-of-mass coordinates: 



N N 
1=1 i=l 



with 



7i 



2n ■ 



(31) 



(32) 



By using Eq. (31) and noting that the last factor of the product depends only on the center- 
of-mass coordinate, which is symmetric under the exchange of nucleons, the A^-nucleon wave 
function can be rewritten as 

^SMsTMri^l: ^n) — ^SMsTMri^l: ^N-i)v's^{'^n) : (33) 

which defines the intrinsic function that depends solely on the relative coordinates, 

*SMsTMt(Si, Sat-i) = A{(fil\{:!Ci), (fisll'^i^N-dXSMs^TMT}- (34) 

As will be shown in the next subsection, this function serves as a generating function of the 
correlated Gaussian basis. 

The Hamiltonian of Eq. (21) can be recast to the relative plus center-of-mass terms as 



N 



H 



'N-l p 2 



P 2 



N 



(35) 



where Pj is the momentum canonically conjugate to the Jacobi coordinate Xj. The matrix 
elements of the Hamiltonian are then 



11 



i^SMsTMri^l, ■-, SAr)|i?|$5'M^T'M^(s'l, s'at)) 

= {^SMsTMri^l: ^N-l)\Hrel\^ S' M'^T' m;^{^' 1: ^V-l) ) (v^s^ Iv^s^,) 

The matrix elements of the Hamiltonian Hrei can be expressed by integrating over Sn 

{^SMsTMri^l: SjV-l) | -f^re/ 1 ^5' M^T'M^, (S'l, S'jv-l)) 

^\2tt) J dSN{^SMsTMT{Sl:-:SN)\H\^s'M'sT'M!^{s'u...,s'N)), (37) 

where use has been made of the single-particle matrix elements in Appendix A and the 
formula 



jdS, (-f (S„ - SV)f e-T(S»-s'„,' ^ (-i)"(2„+ 1)!! g)' . (38) 



Equation (37) shows that the matrix elements between the intrinsic function can be 
obtained in a simple way by factorizing the Sjv-dependent terms from the matrix elements 
in the single-particle basis: The quadratic forms in Eqs. (22), (24) and (25), with the help 
of Eq. (36), should take the form 

sAf = S5f S + 7iv(Siv - SV)' {k = o,p). (39) 

Here the (2iV — 2)-element column vector S is an abbreviation for the set of the Jacobi 
generator-coordinate vectors (Si, Sat.i, S'l, SV-i)- The matrix B^^^ is the (2A^— 2)x 
{2N—2) symmetric matrix defined by dropping the iVth and 2iVth rows and columns of the 
matrix TAf^'^T^ where 

It is clear from Eq. (36) that the polynomials Pi(s, r) and the vectors dj defined in Eq. 

(25) can depend only on the relative generator coordinates. The dependence of the matrix 
elements on the center-of-mass variables, Sat and S'^, can, therefore, be factorized and the 
integration in Eq. (37) reduces to such a simple form as the one in Eq. (38). After the 
integration over Sjv the matrix elements can be expressed in terms of the Jacobi generator 
coordinates in the form: 

1 cr(o)c 



(*5M,TM,(Si, SjV-l)|*5M,TM,(S'i, SVl)) = E ^e'^^^." ^ (41) 

1=1 



N p2 

SMsTMri^li ■■■■> SaT-i)! E TT^ I ^SMsTMt (S'l , S'jv-l)) 

i=l '^f^i 
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N 

(^SMsTMt(Si, Siv-l)| X] ^jil^5'M^T'M^(S'l, SV-l)) 

i<j 

Up 

= / drV{r)e-'''-" Ci^^P.(S, rje-^ssl^^s+D^-r. (43) 

Here the convention of renumbering the set of the vectors (S'l, SV-i) of the ket as 
(Sat, S2N-2) is used to simphfy the notation and thus S stands for the set of the vectors 
(Si, S2Ar-2)- The vector Dj is given by J2'j=i^ D(i-^jSj, where the D{i)j''s are formed from 
the elements of the row vector T>({^T by omitting its A'^th and 2A'^th columns. The column 
vector is defined in Eq. (26). 



C. Integral transformation to the correlated Gaussian basis 

Here we show how to evaluate the matrix elements in the correlated Gaussian basis. Let 
us choose the correlated Gaussians with the form of Eq. (6) and introduce the function 

/xlm(u, X, A) = riKLM{u, x)e-5^^^, (44) 

where 

Af-l 

riKLMiu,x) = v^^+^Y£,m(v), with v = ^ M^Xj = ux, (45) 

i=l 

and where u = (wi, ...,un_i). Note that u is a set of (N — l) real numbers, whereas x are 
the (A^ — 1) Jacobi coordinates. The calculation of the matrix elements becomes simple if 
one uses a generating function of the correlated Gaussian. In fact, the following function g 
is found to be most convenient to generate the function /: 



fKLM{u,:si,A) = - — J dtYiMit) I ^^^^^ g(a,t;u,x, A) j ^ , (46) 

t=\t\=i 

where 

g{a, t; u, x. A) = e-|*^^+"^-t, (47) 

2"n!(2n + 2/ + l)!!- ^ ^ 

Equation (46) is easily proved by using the simple formula 

{y,.t)' = vH' Yl E yumuir. (49) 

n,(>0 m=-l 
2n + l — k 

For a case where the function 6lml{^) of Eq. (7) is needed, the generating function g of Eq. 
(47) must be generalized to include another factor a'V ■ t' . Since the following derivation 
remains essentially unaffected by this generalization, we will assume Eq. (6) as 9lMl{^)- 



L 
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The generating function g can be related to the product of the Gaussians centered around 
{Si;i = 1, A'^— 1} through an integral transformation. To show this, we express the product 
of the Gaussian wave packets as 



N-l 



/detr\ ^ lirx+irSj^-is^rSff 



i=l 



(50) 



with an (N—l) x (A^— 1) diagonal matrix 

/27i ... \ 



272 



V 



(51) 



27iv-i / 



where Sh stands for the set of the generator coordinate vectors (Si, Sat.i). By using the 
familiar formula of the n-dimensional Gaussian integration 



/ 



n\ 3/2 



I det^ / 



(52) 



it is easy to prove the following equation by a direct calculation 

' (detr)3/2 ^ 



g(a,t;u, x. A) 



(47r)(^-i)/Met(r- A)_ 



e iTCT 



(53) 



where t = (Ti, Tjv-i), and 

c = r-^(r-^)r-\ Q = c~^-r, 



(54) 



N-l 

T,^atJ2 (rC)-.^ii,- {i = 1, N-l). (55) 

By combining Eqs. (46) and (53), the correlated Gaussian basis is found to be generated 
from the intrinsic state given in Eq. (34) by the integral transformation 

1 / (detr)-'^/2 \3/2 

A{fKLM{n,^,A)xsMsXTMA = (^(4vr)(^-i)/Met(r-A)J 

X / dtYi^M{i) (^S^^e-^*^^ / dS^e-^^-««-+*«-*5M.TM.(Si, S^_o) ^ • (56) 

t=i 

The matrix elements between the correlated Gaussians are now easily obtained by the 
integral transformation from those expressed in terms of the relative generator coordinates 
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S. Using Eq. (56) gives a general formula to calculate a matrix element for any translation- 
invariant operator O 



{A{fKLM(u, X, A)xSMs^TMT}\0\A'{fK'L'M'(n', X, A')xS'M'gXT'M^}) 

_ 1 / (detr)=^ 

~ BklB'k'l' V(47r)(^-i)det(r- A)det(r- A') 

/. /. ( MK+L+2K'+L' _ 

X / / didi'YMtTY.Mt') [ da^K^Lda,'2K'-,L' ^~''"^'' (57) 

X [dSe-'^^'^^+'^^{^SMsTMA^l,---,SN-l)\0\<i!s'M'T'MdS\,...,^ , , 

^ t=t'=i 

where S = (Si, Sat-i, S'l = Sat, S'jv-i = S2Ar-2) and the matrices, C and Q, and the 
vectors T in Eq. (57), although the same notation is used as in Eqs. (54) and (55), are 
extended to include the corresponding primed quantities of the ket, that is. 




As is shown in Eqs. (41)-(43), the S-dependence of the matrix elements is rather simple 
and the integration over S is done analytically. Since the variables a, a', t, and t' appear 
only through the vector T, those operations with respect to them as implied in Eq. (57) are 
performed systematically. An illustrative example is given in Appendix D. The coupling of 
the orbital and spin angular momenta causes no difficulty. It is very satisfactory aesthetically 
that matrix elements between the basis functions with any sets of the relative coordinates 
can be evaluated in a unified framework without any extra transformation of the coordinates. 
The choice of the set of the relative coordinates amounts to the choice of the matrix A. 



D. Extension to the system of particles of unequal masses 



In this subsection we remark the modifications needed to treat few-particle systems 
containing particles of unequal masses. As you will see, most of the formulation presented 
in subsections 111. A— 111. C remains unchanged. Suppose that the masses of the particles 
are mi, m2, mjv- The width parameter of the Gaussian wave packet is to be changed to 



in 



(59) 



The overlap and Hamiltonian matrix elements are obtained as a function of the generator 
coordinates, (si, Stv, s'^, s'^), in a form similar to the previous case. In this general case 
of unequal masses, the matrix U in Eq. (28) which defines a set of the Jacobi coordinates 
must be generalized to 



/ -1 



U 



mi 
m\2 



mi 



nii2...]V-i 
mi 

mi2...jv 



m2 
mi2 



m.2 





1 



m,i2...]v-i 
mi2...jv 






1 

"'AT 

mi2...iv 



(60) 
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where mi2---i = mi+m2H l-rrii. The reduced masses corresponding to this U are accordingly 

given, instead of Eq. (30), by 

= ^^+l'>^^'l■2■■■^ ^. ^ 1^...^ jv_i^^ and i^n = mu-N- (61) 
mu-i+i 

What is important is to reahze that Eq. (31), most crucial in eliminating the center-of-mass 
motion, still holds even with these modifications as 

N N 

n^::(r.) = nv^^W. (62) 

i=l 1=1 

It is then easy to see that the rest of all the formulae are exactly the same as the case of equal 
masses. We can conclude that the needed modifications noted above are rather trivial and 
simple but still assures the elimination of the center-of-mass motion. As a simple example 
of unequal masses, the system of t + ci + —molecule will be considered in subsection IV.B. 



IV. NUMERICAL RESULTS 

This section is devoted to present the solutions of various N=2 — 7-body problems by 
applying the method described above. To test the method, different potentials (Yukawa, 
Gauss and Coulomb) have been used for bosonic and fermionic systems. Some of the exam- 
ples shown here has its own physical significance, and some other solutions may be considered 
as benchmark test and might be useful in comparison of various few-body methods. One can 
expect that, besides the VMC |^,|] and GFMC [§,|13], other methods will also be extended 
to treat more than = 4-particle systems. As only a few solutions are at present available 
for simple potentials, the examples listed here may help to test other methods. 

As was discussed in subsection II. A, there is no restriction on the choice of the orthogonal 
matrix G. We have found, however, that those special rotation matrices which connect 



different sets of the relative coordinates especially suitable (see also jT^) and will use them 
in what follows. This greatly helps to reduce the number of parameters of G. In the 
stochastic selection of the basis elements, these special matrices G and the parameters of 
the diagonal matrix A' are randomly chosen. The vector u in Eq. (6) is also a variational 
parameter. To avoid an excessively large number of variational parameters we limited u- 
vector values to those which are needed to generate the function 9lml{'^) of Eq. (5) for 
a given set of angular momenta. Comparison of our calculation with others confirms that 
these limitations have not deteriorated the accuracy of the present calculation, that is, our 
trial function is fiexible enough. In the calculation the sets of angular momenta (i.e., the 
sets of these special u vectors) are also randomly chosen. The main advantage of using u 
lies in the simple and systematic evaluation of the matrix elements from the point of view 
of both analytical and numerical calculations. Further test calculations will be needed to 
explore the utility of variational parameter. 

Because the dependence of the matrix elements on the variational parameters is known as 
is shown in Appendix D, one can organize the numerical calculation involved in the random 
search economically. A change of the values of the parameters does not require a recalculation 
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of the whole matrix element. Once they have been calculated for one set of values, to 
calculate them for many more requires virtually no time. The average computational time 
is 10 minutes for a four-body and 2 hours for a six-body calculation on the VPP500 computer 
of RIKEN. 



A. Few-nucleon systems 

We have performed model calculations adopting different central potentials as nucleon- 
nucleon interaction. Some of these model problems have already been solved to high accu- 
racy by various methods and therefore we can directly compare the solutions. The potentials 



used for comparison of different methods are (i) the Malfliet-Tjon (MT-V) potential P9 
which has been most extensively used as benchmark test in few-body calculations, (ii) the 
Volkov "super-soft" core potential ||30|, (iii) the Afnan-Tang S3 (ATS3) potential |^ which 



exhibits a strong repulsive core and incorporates a difference between the spin-singlet and 
spin-triplet channels, (iv) and the Minnesota potential |^ which reproduces the most im- 



portant low-energy nucleon-nucleon phase shifts. The Volkov and MT-V potentials are 
spin-independent, while the ATS3 and Minnesota potentials are spin-dependent. The pa- 
rameters of the interactions are tabulated in Table I. We choose fi^ /m = 41.47 MeV fm^. 
The Coulomb interaction is included only in calculations with the Minnesota potential where 
point charges are assumed and = 1.44 MeV fm. 

The spins (and isospins) of the nucleons are coupled through successive intermediate 
couplings. The spin couplings up to N = 7 nucleons are tabulated in Table II. One nat- 
urally expects and test calculations show that, without spin-isospin coupling, the energy 
convergence is much slower. The number of spin-isospin configurations rapidly increases 
with A^. In the case of ®He, for example, assuming 5" = and T = 1, the wave function 
has 5 X 9 = 45 spin-isospin components. The number of components becomes even higher if 
the interaction has non-central spin-orbit, tensor, etc. parts. The nonlinear parameters are 
not optimized with respect to spin-isospin components, but rather, for each trial choice of 
the matrix A, we select the spin-isospin component that gives the lowest energy. Because 
the matrix elements between different spin-isospin components differ only in linear factors 
(cre in Eq. (20)), the calculation of the matrix elements of each spin-isospin component of 
the wave function requires essentially the same computational effort as that needed by the 
calculation for only one component. 

Each calculation has been repeated several times starting from different random points 
to check the energy convergence. The energy as a function of the number of basis states 
is shown in Fig. 2 for the case of ^Li with the Volkov potential. The energies on different 
random paths, after a few initial steps, approach to each other and converge to the final 
solution. The energy difference between two random paths as well as the tangent of the 
curves give us some information on the accuracy of the method on a given size of the basis. 
The root mean square (rms) radius of the few-nucleon system is calculated in each step and 
found to be rapidly convergent to its final value. By increasing the basis size the results can 
be arbitrarily improved when needed. 

The number of basis states required to reach energy convergence increases with the 
number of particles but it depends on the form of the interaction as well. This latter 
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property is illustrated in Fig. 3 for the case of the a-particle. The soft-core Volkov potential 
shows rapid convergence, while the hard-core ATS3 interaction requires more basis states 
to get an accurate solution. The relatively fast convergence for the MT-V potential of a 
strong repulsion can be explained by the simplicity - the spin-independent nature - of this 
interaction. 

In the following we show tables for the ground-state energies E and point matter rms radii 
(r^)^/^. The basis dimension /C of the SVM listed in the tables is such that, beyond it, the 
energy and the radius do not change in the digits shown. Table III shows our results (SVM), 
together with results of others, for the application of the spin-averaged MT-V potential 
to A^=2-7-nucleon systems. For three-body systems, the solution of the Faddeev equation is 
known to be the method of choice, but the SVM can easily yield the same accurate energy. 
As the MT-V potential is a preferred benchmark test of the few-body calculations, there 
are numerous solutions available. Table III includes a few of the most accurate results. 
The nice agreement for four-nucleon case corroborates that the SVM is as accurate as the 
direct solution of the FY equations |33], the method of the Amalgamation of Two-body 



correlations into Muhiple Scattering (ATMS) @] process or the VMC 0] and GFMC 
method. The basis used in the Coupled-Rearrangement-Channel Gaussian-basis (CRCG) 
variational method |^ is similar to that of the SVM, but the Gaussian parameters follow 
geometric progressions. The fact that the basis size needed in the SVM is much smaller 
proves the efficiency of the selection procedure. The results of the VMC calculation for the 
five- and six-nucleon systems are also in good agreement with the results of SVM. The MT- 
V potential has no exchange term; therefore, unlike the nature, it renders the five-nucleon 
system bound, and the nucleus tends to collapse as the binding energy increases with the 
number of particles. 

The next example is the Volkov potential which, due to its very soft core, is the most 
readily solvable case. This simple potential is widely used in model calculations for light 
nuclei. As one sees in Table IV, the results of SVM agree with those of other calculations, 
especially with the one using hyperspherical harmonics (HH) functions. The number of 
basis states needed to reach convergence is remarkably smaller than in the case of the MT-V 
potential. Without the Majorana exchange term (M = 0) this potential also leads to a 
collapsing system. By setting the Majorana parameter to M=0.6, a commonly used value 
to get the correct binding energy, one may obtain more reasonable energies. The Volkov 
potential with M=0.6 does not change the energies of N=2—4 nuclei, does not bind ^Li in 
accordance with the nature, but does bind the ^Li ground state {E=— 31.82 MeV, (r^)^^'^ 
=2.69 fm). 

Another potential that is often used in test calculation is the ATS3 potential. We have 
challenged the SVM to get solution for this case because, unlike the Volkov, this spin- 
dependent potential has a relatively strong repulsive core (see Table I). The solution, al- 
though on a somewhat larger basis, can easily be obtained, and it is in good agreement 
with those of other methods in the N=3 system as shown in Table V. We note, however, 
that the energy of the SVM is significantly lower than the ones of other methods for the 
a-particle. This may be due to the strong repulsion of this potential. For example, the FY 
calculation [^] agrees with ours for the MT-V potential, but shows a noticeable difference 
in the case of the ATS3 potential. The variational calculation |jT6| using a correlation factor 
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also misses considerable energy for the ATS3 potential. Surprisingly, inspite of its exchange 
part, this potential also binds the five-nucleon system and overbinds the six-nucleon systems 
very much. 

The last example for the few-nucleon system uses the more realistic Minnesota poten- 



tial |22[, which is a central interaction of Gaussian form, containing space-, spin-, and 
isosp in-exchange operators (see Table I). The Minnesota potential has often been used in 
cluster-model calculations of light nuclei. Table VI shows results with this potential, where 
the Coulomb interaction between protons is also included. All possible spin and isospin 
configurations are allowed for and all spherical harmonics that give non-negligible contri- 
bution are included in 6'lj\/^(x). Since the method has proved to be accurate and reliable 
for other potentials, it is justifiable to view these results as testing the interaction rather 
than the method, and hence the results are compared with experimental data. The energy 
and the radius of triton and a-particle converge, with small bases, to realistic values. The 
Minnesota potential, correctly, does not bind the A^=5 system, but it binds ^He and slightly 
overbinds ^Li. The radius of ^He is found to be much larger than that of ^He, consistently 
with the halo structure of ^He . It is for the first time that the Minnesota force is tested 
without assuming any cluster structure or restricting the model space by any other bias. 
The agreement is surprisingly good not only with experiment but also with cluster-model 
calculations for all nuclei [p^] . 

It is interesting to note that none of these simple potentials binds the four-neutron 
system. The Volkov potential, for example, is so strong that it binds the singlet two-neutron 
system, but it does not allow the neutrons to form a four-neutron bound state due to the 
Pauli principle. 

As an example for bosonic nuclear few-body system, we consider the case of structureless 
a-particles interacting via the state-independent potential; 

V{r) = 500 exp [-{0.7rf ] - 130 exp [-{0A75rf ] (MeV), (63) 



where r is in fm. This potential is taken from Ali-Bodmer's S'-state potential It has 



a repulsive core which is about 370 MeV high and extends up to 2 fm. The repulsive 
core prevents the a-particles from collapsing. The results are compared in Table VII. Our 
calculation agrees with the ATMS result for the N =3 and 4 systems. 



B. Coulombic systems 

The results for the long range 1/r potential are collected in this subsection. The first 
example is the polyelectric system {■me~^ ,ne~). The possibility that m positrons and n 
electrons form a bound system was originally suggested by Wheeler and this question has 
been extensively studied since then. Besides the trivial and analytically solvable m = 1 and 
n = 1 case, the existence of the positronium negative ion (m = 1, = 2) was also predicted 
by Wheeler [|D|. Dozens of works have attempted to solve the e"'",e~,e~ Coulombic three- 



body problem, continuously refining the accuracy of the calculated binding energy pT| , ^ 



Despite of numerous attempts, no one has obtained bound states for the polyelectric system 
of more than four particles. The positronium negative ion has experimentally been observed 
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The binding energy of the positronium molecule (2e"'", 2e ) was first calculated by Hyller- 



aas and Ore |^4[. To date, the positronium molecule has not been directly observed, and 
this fact intensifies the theoretical interest to solve this Coulombic four-body problem. This 
molecule is short-lived because the electron and positron may annihilate. Unlike the positro- 
nium ion, the positronium molecule is neutral, and therefore the best chance to distinguish 
it from the positronium itself is related to their different lifetimes. The QED formulae 
to determine the probability of a pair annihilation in the positronium molecule through a 
/c-photon process {k = 0,1,2,...) would require a highly accurate wave function [^ . 

In Table VIII we compare our results to the most precise calculations found in the 
literature. The correlated Gaussian function without the polynomial part {K = in Eq. 
(6)) is known to poorly represent the Coulomb cusps HTSf. To improve the cusp properties 
the trial function with = 0, 1, 2 polynomials has been used. 

As is shown in Table VIII, our calculation reproduces the first six digits of the variational 
calculation of Ref. |5T] for the ground state of (2e"'",e~), and the rms radius also agrees with 
it. There are two recent variational calculations [[19|,^ for the positronium molecule using 
the correlated Gaussian functions. In these works the nonlinear parameters were determined 
by optimatization. To compare our calculations to theirs directly, the value of K = was 
chosen and the same basis size (/C=300) was used. Our result is slightly better than the 
energy obtained by them and this reinforces the reliability and powerfulness of the random 
selection of the nonlinear parameters. The number of nonlinear parameters of this case is 
/Cx(4x3)/2 = 1800. The complete optimatization of the parameters is, of course, superior 
to the SVM. Test calculations show that, provided the number of parameters is low, that 
is, a full optimatization is feasible, the optimatization finds lower energy. But when the 
number of parameters becomes high, the full optimatization becomes less and less practical 
partly because it fails to find true minimum and partly because the computation becomes 
too excessive. 

We found no bound states for the (3e"^,2e~) and (3e"'",3e~) systems. The energy of 
(3e"'",3e~), for example, converges to the sum of the energy of a dipositronium molecule 
and of a positronium (0.515989 a.u.+0.25 a.u. =0.765989 a.u.). Allowing the selection of the 
parameters from a larger region increases the rms radius, which is typical of an unbound 
state. The system of a negative and a positive positronium ion thus forms no bound state 
but dissociates into a dipositronium molecule and a positronium. 

Calculations for L = 1 state also fails to find a bound system. This result entails that 
the Coulomb force cannot bind more than four particles out of identical charged fermions 
and their antiparticles. 

To examine the role of the Pauli principle in preventing five-electron-positron system from 
forming bound states, we repeated the same calculation replacing the fermions by bosonic 
equivalents. On a different scale, these systems may be identified, e.g., as the systems of 7r~ 
and vr"*" with their strong interaction neglected |^^. Such bosons turn out to form bound 
states even for = 5. As is expected, the radius of the charged boson system decreases 
as the number of particles increases. It is interesting to note that the energies of bosonic 
and fermionic systems are equal for A^ = 3 and for A^ = 4. The reason is that the energy 
minimum belongs to the same spatial configurations, that is, to a triangular pyramid for 
A^ = 4, for example |^ . 



In Table IX we show results for bosonic and fermionic systems with a purely attractive 
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Gm?/r ("gravitational") interaction. An A^-body system of identical particles bound to- 
gether by attractive pair potentials always collapses in large-iV limit (the binding energy 
per particle rises with N to infinity), even if the particles are fermions. Self-gravitating 
boson systems have recently attracted some interest For these systems, both varia- 

tional lower and upper bounds are available. In this case even the five-fermion system is 
bound. Thus the lack of bound states in five-electron-positron systems is a joint effect of 
the antisymmetry and of the repulsion between identical particles. 

Finally, we mention an example involving an excited state. With K, = 500, the SVM 
gives the energies of the ground and first excited states of the t + d + fi~ system as —111.3640 
and —100.9121 a.u., which are respectively compared to —111.364342 and —100.916421 a.u. 
of the CRCG result with JC = 1442, while the configuration-space Faddeev calculation 
H] gives —111.36 a.u. for the ground state but no information for the excited state. 



V. SUMMARY 

We have formulated a variational calculation for few-body systems using the stochastic 
variational method on the correlated Gaussian basis. We have demonstrated the versatility 
of the correlated Gaussians and the efficiency of the stochastic variation by various numerical 
examples for A^ = 2 — 7-particle systems. All the details of both formulation and calculational 
procedure are included to make this paper self-contained and easily reproducible. 

The comparison with other calculations has corroborated the accuracy and efficiency of 
the method. In none of the test cases has the present method proved to be inferior to any 
of the alternative methods, and yet the method does not require excessive computational 
effort. 

The correlation between the particles plays an important role in describing the few- 
body system realistically. It has been taken into account in the framework of the correlated 
Gaussian functions. The correlated Gaussians are constructed from products of the Gaussian 
wave-packet single-particle functions through an integral transformation, which has enabled 
us to evaluate the center-of-mass motion free matrix elements analytically starting from 
the single-particle level. The nonlinear parameters of the correlated Gaussians have been 
selected by the stochastic variational method with a trial and error procedure. The success 
of the method using the correlated Gaussian basis is probably due to the fact that none of 
the Gaussians is indispensable, that is, there are different sets of the Gaussian parameters 
that represent the wave function equally well. 

The method presented in this paper can be useful to solve few-body problems in diverse 
fields of physics such as description of microclusters, non-relativistic quark model, and halo 
nuclei. Among others, the most important application is the solution of the nuclear few-body 
problem, that is, a description of light nuclei by using realistic nucleon-nucleon potentials. 
In this case one has to take into account both short-range repulsion and higher orbital 
angular momenta required by the non-central components. Our test examples show that 
the correlated Gaussian basis function might be a suitable candidate to cope with these 
requirements. As is explained in Appendix B, the evaluation of the matrix elements for the 
non-central potentials poses no serious problem, and calculations including such potentials 
for nuclear few-body systems are under way. 
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The limitations of the present method are those imphed by the basis size and by the 
computer memory to store the matrix elements in the generator coordinate space. The 
limitations may become excessive as the number of particles and/or spin and isospin con- 
figurations become large. 

To extend the method to nuclei of larger mass number in an approximate way, one can 
freeze part of the model space for a group of nucleons (cluster). One can omit, for example, 
some of the spin or isospin channels. It might be a good approximation to consider only 
those spin channels where the spins of the like nucleon pairs are coupled to zero. One can 
also restrict the intrinsic spatial motion of a cluster by fixing the nonlinear parameters to 
some appropriate values. One can introduce N clusters and place the nucleons of each cluster 
into a common harmonic oscillator well, for example. The microscopic multicluster model 
is based on this approximation. The matrix elements needed in this multicluster model are 
given as a special case of those presented in this paper. In fact, one only needs to choose 
the single-particle generator coordinates, (si, ...,8^), such that 

Sn,_i+1 = S„,_i+2 = ■ ■ ■ = s„^_^+„^ (i = 1, AT), (64) 

where is the number of nucleons in the ith cluster (no = 0, rii + ■ ■ ■ + un = A), and 
needs to couple appropriately the spins and isospins of the nucleons in each cluster to the 
spin and isospin of the cluster. The microscopic multicluster model has been successfully 
applied for description of the structure of light nuclei (see, for example, 1 17,22,23]). 



Finally, we summarize some merits of our method in the following: 

(i) Fully analytical calculational scheme; this plays a major role in the high speed and 

accuracy of the calculation. 

(ii) Universality of the scheme. One needs to introduce no change, for example, between 
describing a multinucleon system and a Coulombic few-body system. It is easily adapt- 
able to identical or non-identical particles, to fermions or bosons or mixed systems. 
The masses of particles may be different, yet no problem with the center-of-mass mo- 
tion arises. 

(iii) No expansion of the interaction is needed, and thus no problems in partial-wave trun- 
cation arise. 

(iv) The convergence of the energy is fast. If one needs just a 2-3-digit estimate of the 
energy, it is enough to use a very small basis. 

(v) The method is also accurate for excited states, which are obtained simultaneously with 

the same diagonalization (provided their angular momenta and parities are the same 
as those of the ground state; but only such excited states may be problematic). 

(vi) The wave function is obtained in a compact, analytical form. It is then easy to use it in 
calculations of physical properties. It is "portable" , reproducible, and easily testable. 
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APPENDIX A. SINGLE-PARTICLE MATRIX ELEMENTS 



The aim of this appendix is to hst the single-particle matrix elements between Gaussian 
wave packets (Eq. (17)). The Gaussian packets are generalized in this appendix to have 
different width parameters and the expressions are therefore slightly more general than 
needed in the formulae of the main text. These single-particle matrix elements are, however, 
required for treating particles of unequal masses, in which the width parameter v belonging 
to the particle of mass m is to be chosen by Eq. (59). The overlap of two Gaussian wave 
packets is 




The matrix element of the kinetic energy operator (T = — ^A) reads as 



The matrix element of the square radius becomes 

1 / 2 



The two-particle matrix element of a 5-function is given by 



(«|5(ri-r2-r)|«) 



^(z/i + ^^3 + '^4)7r / 

X exp ■ (r ■ + ■ j W^s ¥^S3 W^S2 y^sl • 68 

\ Z/l -I- i/2 + + ^ Vl^VZ Z/2 Z/4 ^ / 



APPENDIX B. MATRIX ELEMENTS OF THE TWO-BODY POTENTIALS 

The scope of this appendix is the calculation of the matrix elements of the different 
ingredients of the two-nucleon interaction between Gaussian wave packets. Most of the 
widely used coordinate-space two-nucleon interactions consist of central- (OJ2), tensor- (OJ2), 
spm- orbit- (OJ2), (or p^) (O12) and quadratic spin-orbit (OJ2) type potentials. (We follow 
the abbreviated notations - c,t,b,q,bb - invented by Urbana-Argonne group The 
definition of these operators is given in Table X.) These potential terms might be multiplied 
by the Ti-T2 {OI2), cri-cr2 or Ti-T2 cri (T2 {OI2) spin-isospin operators and then one 

ends up with the general form 

Vi, = J2[ dvVPir)6ir, - r, - r)Of^, (69) 

p 

where p = c, ct, ca, era, t, tr, b, br, q, qr, qa, qra, bb, bbr is the short-hand notation to specify 
the component and V^{r) is the corresponding radial form factor. By using Eqs. (16), 
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(17), (68), (69), and after some straightforward transformation, the matrix element of this 
interaction can be written as 

(C..nC..^l^i2lCv;r(C'.^.^) = E / drV^{r)f,{r)M{v) 

p •' 

X X^,.,X.Jr'B^ + E E (-l)Vl^^(r)Cf_^ | x^.^X.^^ ^H^^)^ (70) 

1=0 m=-l 

with 

(71) 

The and Cf^, independent of r, are the operators in spin-isospin space and are hsted 
in Table X for the most important terms {c,t,b,q,bb). The remaining terms can be easily 
derived by multiplying these operators by the appropriate <t and r operators. The function 
fp{r) has the simple form: fp{r) — ii p — t or tr, and fp{r) — 1 otherwise. 

The calculation of matrix elements of the operators appearing in Table X in the spin part 
is easily done with the Clebsch-Gordan coefficient. They are given below by suppressing the 
spin function x^ai i^i = 

(mi 7712 I O"! • <7"2 I rn'i ) 

I 1 11 1 

= 3(-l)™l""''l5^,+m2X,+m' I T: mi) {-1712 I 1712-^2 I o"^2)- (72) 

(mi m2 I [cTi X (T2]J^^ I m[ m'a ) = 3 (5m,mi+m2-ml-mi 

X (- m'l 1 mi^m'i | - mi) (- m2 1 m2—m'2 \ - m2) (1 mi^m'i 1 m2-^2 | 2 m) . (73) 



(mi m2 I [x X (cTi + CTs)]^^ | m[ m'2) ^ (x)^^ 

'?l,'?2 = -l 

X (1 5i 1 52 I 1 "^){5m2,m^(^ 1 ?2 I ^ mi) + (5mi,m'i(^ "^2 ^ ^2 | ^ ^2)}. (74) 

(mi m2 I (x • o-2)<ri^ + (x • (Ti)<r2^ | m[ m'2 ) 

= 3 {(-1)-^— 2 5m,m,-m',{^)m',-m, + {-1)^'-^'^ 5m,m,-m',{^)m',-m,} (75) 

X (- m'l 1 mi— m'l | - mi) (- m'2 1 1712-^2 \ - ^2). 



(mi m2 I (x • <Ti)(x • 0-2) I m'l m'2) = 3 (_i)'»i-S+'»2-mi 
X (x)^'^_^^(x)„^_„2(- m'l 1 mi-m'i | - mi){- m'2 1 m2^'2 | - ^2). (76) 
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(mi 777,2 I [x X <Ti]^^^ X [x X <T2]^^^ I 777 i 7772 ) 

= 3 (- m[ 1 TTii^m'^ I - TTii) (- Trig 1 1712-^2 \ - ^2) 
1 

^ H (x)gi(x)g2(lgi+mi-n7'i lg2+?«2^2 1 2 m) (77) 

qi,«2=-l 

X (1 Q'l 1 mi^m'i I 1 5i+mi^mi) (1 Q'2 1 m2^/772 | 1 g2+n72^?772) . 
Here x is a 3- dimensional vector and {'x)m stands for its spherical component ^J^xYijn{x). 

APPENDIX C. MATRIX ELEMENTS OF SLATER DETERMINANTS 

In this appendix we briefly outline the calculation of the matrix elements between Slater 
determinants of Eq. (18) and show their concrete functional form, that is, the dependence 
on the generator coordinates s. We assume that the width parameter of the Gaussian wave 
packet is chosen to be the same for all nucleons. The overlap of two Slater determinants is 
equal to the determinant of the matrix of the single-particle overlaps: 

(0,(si,...,Sjv)|0«'(s'i,-,sW)) = det{5}, (78) 

where 

B,j = {z,j = 1, N). (79) 

J J 

By using the deflnition of the determinant, this can be rewritten as 

iV! 

{M^i,...,s^)\M^\,...,s'^))^Y.^ign{Pm^^^^^^ (80) 

where {pi ■ ■ ■ Pn) is the permutation P of the set (1 ■ ■ ■ N). Substituting the overlap of the 
single-particle overlaps of Eq. (65) into Eq. (80) yields an explicit formula for the overlap 
of the Slater determinants: 

AT! 

((/)„(si,...,s^r)|0«'(s'i,...,s;,)) = 5:Cpe-^^^-^ (81) 

p 

where the matrix is deflned by 

{Ap)ij = {Ap)N+i,N+j = i^Sij, {Ap)i^N+j = {Ap)N+j,i = -i^Sj^p^ (1 <i,j <N), (82) 
and 

Cp = Sign(P)5^,^,^^ 6r,r'^, ■ ■ ■ ^a^a'^^ ^r^r'^^ ■ (83) 

The orthogonality of the spin-isospin functions greatly reduces the number of terms in the 
summation over P. 

The matrix element of the kinetic energy operator is 
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N N N 

EiM^l^---,^N)mA^\,...,s'^))=Y:Ei^s.a^^^^^^ (84) 

1 • 1 • 1 ■'3 1 

1=1 1=1 ;?=1 

where B^^ is obtained by omitting the ith row and the jth column of the matrix B defined 
in Eq. (79). The substitution of the single-particle overlaps and the single-particle matrix 
elements of the kinetic energy operator enables us to obtain 

N hi ) ^- /"^ ^ \ - 

x:(0«(si, s^m\M^\, s'^)) = ^ E (i^ - 2^^^^) (85) 

Note that the coefficients Cp and the matrices Ap in Eqs. (81) and (85) are the same. The 
manipulation similar to this leads us to the matrix element of the square radius as 

E(0.(si, s^,)|r2|0,,(s^, s'^)) = -Y,Cp(^-N- -sAps + u~ssj e-^^^-^ (86) 

As is explained in Appendix B, the matrix element of any two-body interaction between 
the Slater determinants may be reduced to the following: 

f;(0,(si, ...,Sr,)\S{r, - rj - r)Ofj\Ms\, ...,s'^)). (87) 

i<j 

The matrix element of Eq. (87) is called the correlation function of type p evaluated between 
the Slater determinants. The calculation of this matrix element can be done with the use of 
the basic two-body matrix elements of Appendix B and the single-particle overlaps. We will 
show, as an example, the case of Ofj = 1, i.e., Wigner-type 5-function two-body interaction 
{p — c). Then we have 

N 

E(0«(si, sjv)|5(ri - Tj - r)|(/)«/(s'i, s'^)) 

i<j 

N N 

= E E {0s...rM..,rJS{r. - - r)|^.\.,..^r,.;.;)(-l)^+^-+'=+'det{i?^^'='}, (88) 

i<j k,l=l 

where S*^*^' is obtained by omitting the ith and jth rows and the kth and Zth columns of 
the matrix B. By substituting the explicit formula of the ingredients, we arrive at the same 
form as given in Eq. (25): 

N 

E(0«(si, SN)\S{ri - Tj - r)|(/)«/(s'i, s'jv)) 

i<j 

N< N 

^ (^)3/2e— ^ E E e-i^(^-+^^"')^+^^^'-, (89) 

with 

dj^-?) = u{Si - Sj + SN+pi - Sn+pj). (90) 

26 



The coefficients Cp and the matrices Ap are the same as those that appear in Eq. (81). 
The matrix Sp"''* is defined by 



K^P Jkl \^P Ilk \ _^ — (A M_L^.\ (n M^^.\ (M^^. AT^^.\ V»-^y 



2 (kl) = {tj),{t,N + pj),{j,N + pi),{N + pi,N + pj). 
otherwise. 



APPENDIX D. EVALUATION OF THE OVERLAP OF THE CORRELATED 

GAUSSIANS 

In this appendix we derive the overlap of the basis functions, defined in Eqs. (2) and 

(6), to illustrate the calculation of the matrix elements through the operations prescribed 
in Eq. (57). For the sake of simplicity we work out the case of overlap only, but the matrix 
elements of the kinetic and the potential energy are not much more involved either and these 
matrix elements can be calculated by repeating the steps detailed here. 

In Eq. (57) the matrix element of an operator O between the correlated Gaussians are 
derived from the matrix elements between the Gaussian packets by applying suitably chosen 
operations. In the case of (9 = 1, upon substituting the matrix elements in the generator 
coordinate space (Eq. (41)), Eq. (57) reads as 

(^{/kLm(u, X, A)xSMs^TMT}\A'{fK'L'M'{u', X, A')xS'M'gXT'M^}) 

11 didi'YLMiiyYL'M'ii') (92) 

/ JK+k' ^ . "o ^ 

\ da'^da"^ J 



i=l / "="'=0 

t=t'=i 



with 



K = 2K + L, k' = 2K' + L', (93) 

where a constant factor is omitted for a moment. By integrating over the generator coordi- 
nates S with the use of Eq. (52), one readily obtains 

(^{/kLm(u, X, A)xSMsXTMT}\^{fK'UM'W, X, >1')X5'M^'^T'M^}) 

didi'YLM{irYL,M'{t') (94) 



X 



d-+-' ^^(o)f (27r)(2^-2) 



da'^da'-' ^ ' \det{B\"^ +Q) 



a=cx =0 
t=t'=l 



where, by using Eqs. (54) and (55), we introduced the abbreviations 

1 7V-1 1 27V-2 

= ^ E ^AiBt^ + Q)-' - c)jkWk, P/ = ^ E + Q)'' - c)jkWk', 

j,k=l ^ j,k=N 

N-1 2N-2 

= E E ^jUB^"^ + Q)~' - C)jkW, (95) 

j=l k=N 



27 



with 



N-l 2N-2 

^fe = E (rC)^^«„ w,' = Y: i^C%lu/, (96) 

j=l j=N 

in order to emphasize the a- and a'-dcpcndcncc of the resulting expression. The differentia- 
tion with respect to a and a' can be promptly given by expanding the exponential functions 
into power series: 

^K+K ^ ^^2_^_p.f^,2_^_g.^^,^,^f 



da'^da"^ 

j\j'\k\ (2j + k-K)\i2j' + k-K')\ ' ^ ' 

and therefore, after putting a = a' = 0, the calculation of Eq. (94) can be continued as 

- 11 didi'YLMiiyYL'M'ii') (98) 



X 



t=t' = l ' 



where the summation over k runs from to max(fi;, k.') for those values that fulfill (—1)*^+'^ = 1 
and (-1)'=+'^' = 1. 

The last step, the integration over the angles of t and t', can be accomplished by applying 
Eq. (49) to express the scalar product t • t', and then we get the final expression: 

{A{fKLM{u, X, A)xSMs^TMT}\A'{fK'L'M'{u', X, A')X5'M^'Vr'M^}) 

1 / (detr)3 
^ '^""^^^'S^.fi^,,, i(4v^)(^^-l)det(^-A)det(^-A') j ^^^^ 



X 



where the constant term previously suppressed is also included. The orthogonality of the 
matrix element in the spin and isospin quantum numbers is implicitly contained in the 

coefficients C^"^ . 

The calculation of the matrix elements in the way described above is very simple. This 
fact becomes especially evident if one compares the work described above to the formidable 
task of the calculation of the matrix elements in the case where the function 9lml{'^) is 
decomposed into partial waves of the relative coordinates as in Eq. (5). In fact, in that 
latter case one has to integrate over the angles of the relative coordinates and one has to 
cope with complicated angular momentum algebra. We note, however, that the calculation 
of the matrix element of the latter type poses no problem if the function 9lml{'^) of Eq. (5) 
is expressed as a linear combination of the terms of Eq. (6) with appropriate u- vectors. 
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All the matrix elements can be given in a similar closed analytic form and the numerical 
evaluation of the matrix elements as a function of the nonlinear parameters is therefore 
straightforward. The values k = 2K+L and k' = 2K'-\-L' are usually small in practical cases 
and the summation over k is limited to just a few terms in Eq. (99). 
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TABLES 



TABLE L List of the parameters of the nucleon-nucleon potentials used in this paper. The 
potential consists of a few terms; each is expressed as Vif{fii,r){Wi + BiP^j — HiP-^ + MiPr), where 
Po-, Pr and Pr are the spin-, isospin- and space-exchange operators. The form factor f{fii,r) is 
exp(— ^r^) for the Gaussian potential or exp(— //r)/r for the Yukawa potential. The potential 
strength V is in units of MeV and the range /x in units of fm~^ for Gaussian or fm~^ for Yukawa, 
respectively. The Majorana mixture M of the Volokov potential is set to zero in the calculation. 
The Minnesota potential contains the parameter n, which is set to unity in the calculation. 



Potential 


Type 


i 


Vi 






Mi 


Bi 




MT-V 


Yukawa 


1 


1458.05 


3.11 


1.0 


0.0 


0.0 


0.0 


m 




2 


-578.09 


1.55 


1.0 


0.0 


0.0 


0.0 


Volkov 


Gauss 


1 


144.86 


0.82-2 


1.0 -M 


M 


0.0 


0.0 


m 




2 


-83.34 


1.60-2 


1.0- M 


M 


0.0 


0.0 


ATS3 


Gauss 


1 


1000.0 


3.0 


1.0 


0.0 


0.0 


0.0 


m 




2 


-326.7 


1.05 


0.5 


0.0 


0.5 


0.0 






3 


-166.0 


0.80 


0.5 


0.0 


-0.5 


0.0 






4 


-43.0 


0.60 


0.5 


0.0 


0.5 


0.0 






5 


-23.0 


0.40 


0.5 


0.0 


-0.5 


0.0 


Minnesota 


Gauss 


1 


200.0 


1.487 


0.5u 


1.0-0.5ti 


0.0 


0.0 


m 




2 


-178.0 


0.639 


0.25^ 


0.5-0.25ii 


0.25n 


0.5-0.25^ 






3 


-91.85 


0.465 


0.25^ 


0.5-0.25^ 


-0.25n 


-0.5 + 0.25^ 
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TABLE II. Spin (isospin) couplings of the micleons. The possible intermediate spin (isospin) 
configurations of an A^-nucleon system are listed above the (A^— l)th horizontal line. The spins of 
a 5-nucleon system, for example, are given by the first ten elements of the 5th column. The spins 
of the 2 — 3 — 4- particle subsystems of the 5-nucleon system is given in the first ten elements of 
the 1 — 2 — 3 columns. 



N 


2 


3 


4 


5 


6 


7 







1/2 





1/2 





1/2 




1 


1/2 





1/2 





1/2 




1 


3/2 


1 


1/2 





1/2 







1/2 


1 


1/2 





1/2 




1 


1/2 


1 


1/2 





1/2 




1 


3/2 


2 


3/2 




1/2 




1 


3/2 


1 


3/2 




1/2 







1/2 


1 


3/2 




1/2 




1 


1/2 


1 


3/2 




1/2 




1 


3/2 


2 


■j/2 




3/2 







1/2 





1/2 




1/2 




1 


1/2 





1/2 




1/2 




1 


3/2 


1 


1/2 




1/2 







1/2 


1 


1/2 




1/2 




1 


1/2 


1 


1/2 




1/2 




1 


3/2 


2 


3/2 


2 


3/2 




1 


3/2 


1 


3/2 


2 


3/2 







1/2 


1 


3/2 


2 


3/2 




1 


1/2 


1 


3/2 


2 


3/2 




1 


3/2 


2 


5/2 


3 


5/2 




1 


3/2 


2 


3/2 


1 


3/2 




1 


3/2 


1 


3/2 


1 


3/2 







1/2 


1 


3/2 


1 


3/2 




1 


1/2 


1 


3/2 


1 


3/2 




1 


3/2 


2 


5/2 


2 


5/3 







1/2 





1/2 


1 


3/2 




1 


1/2 





1/2 


1 


3/2 




1 


3/2 


1 


1/2 


1 


3/2 







1/2 


1 


1/2 


1 


3/2 




1 


1/2 


1 


1/2 


1 


3/2 




1 


3/2 


2 


3/2 


2 


5/2 




1 


3/2 


1 


3/2 


2 


5/2 







1/2 


1 


3/2 


2 


5/2 




1 


1/2 


1 


3/2 


2 


5/2 




1 


3/2 


2 


5/2 


3 


7/2 
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TABLE III. Energies and rms radii of A'^-nucleon systems interacting via the Malfliet-Tjon po- 



tential V [29|. The value of K, in Tables III— IX denotes the basis dimension beyond which the en- 
ergies and the radii of the SVM calculation do not change in the digits shown, {fi^ /m = 41.47 MeV 
fm^ throughout Tables III— VI). 



N 


{L,S)J^ 


Method 


E (MeV) 


(r2)i/2 (fm) 


K 


2 


(0,1)1+ 


Numerical 


-0.4107 


3.743 








SVM 


-0.4107 


3.743 


5 



(0,l/2)l/2^ 



Faddeev ^ 
ATMS II 
CHH g 
GFMC [i 



VMC II'' 
SVM 



.25273 

.26 ±0.01 

.240 

.26 ±0.01 
.2689 ±0.03 
.2527 



1.682 

1.682 

1.68 

1.682 



80 



(0,0)0^ 



FY H 
ATMS II 
CRCG § 
GFMC [i 



VMC 
SVM 



-31.36 

-31.36 

-31.357 

-31.3 

-31.3 

-31.360 



±0.2 
±0.05 



1.40 

1.36 
1.39 
1.4087 



1000 



150 



(1,1/2)3/2- 



VMC |34| 
SVM 



-42.98 
-43.48 



±0.16 



1.51 
1.51 



500 



6 C'He) 



(0,0)0^ 



VMC 

SVM 



-66.34 
-66.30 



±0.29 



1.50 
1.52 



800 



7 (^Li) 



1/2)3/2- 



SVM 



-83.4 



1.68 



1300 



^Calculated with Coulomb potential, the Coulomb contribution then subtracted perturbatively. 
The potential strength used in the VMC @ calculation is ^1=1458. 25 and V2 = -578.17 MeV, 
which is slightly different from that used in the present calculation. 
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TABLE IV. Energies and rms radii of A'^-nucleon systems interacting 
|30|. The Majorana exchange parameter M is set to zero. 


via the Volkov potential 


N 


{L,S)J^ 


Method 


E (MeV) 


(r2)V2 (fm) 


K 


2 


(0,1)1+ 


Numerical 
SVM 


-0.545 
-0.545 


3.44 
3.44 


5 


3 


(0,1/2)1/2+ 


Fadd 
Varic 
HH 1 
SVM 


eev II 
itional ||l^ 
35| 


-8.43 
-8.460 
-8.4647 
-8.46 


1.725 
1.73 


30 


4 


{ c\ f\\f\-\- 
(0,0)0^ 


FY J 

Varic 
HH 1 
SVM 


1 

itional [l^ 
36 1 


— 6\J.I( 

-29.490 

-30.3988 

-30.42 


1.47 
1.49 


50 


5 


(1,1/2)3/2- 


SVM 


-43.00 


1.59 


120 


6 (6Li) 


(0,1)1+ 


SVM 


-66.25 


1 60 


250 


7 (^Li) 


(1,1/2)3/2- 


SVM 


-98.75 


1.57 


400 


TABLE V. Energies and 
potential plj. 


rms radii of A^-nucleon systems interacting 


; via the Afnan-Tanj 


I S3 


N 




Method 


E (MeV) 


(r2)V2 (fm) 


/C 


2 


(0,1)1+ 


Numerical 
SVM 


-2.216 
-2.216 


1.94 
1.94 


7 


3 


(0,1/2)1/2+ 


Faddeev Q 
Faddeev @ 
GFMC [10| 
Variational [[l6[ 
SVM 


-8.20 
-8.765 

-8.73 ±0.10 

-6.677 

-8.753 


1.72 

1 797 

1.67 


40 


4 


(0,0)0+ 


FY II 

Variational 
SVM 


-28.80 

-25.654 

-30.37 


1.44 
1.42 


140 


5 


(1,1/2)3/2- 


SVM 




-44.27 


1.48 


220 


6(6Li) 


(0,1)1+ 


SVM 




-70.65 


1.49 


400 
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TABLE VI. Energies and rms radii of A'^-nucleon systems interacting via the Minnesota po- 
tential [^] with the exchange parameter u = 1. The Coulomb interaction is included (e^ = 1.44 
MeV fm). The experimental (r^)^/^ value is the point charge radius with the proton's finite size 
corrected. 



AT 




iVitJLllUQ 


H/ I ivie V ) 


(r2)i/2 (frn) K. 


2 


(0,1)1+ 


SVM 
Exp. 


-2.202 
-2.224 


1.952 
1.96 


5 


3(3h) 


(0,1/2)1/2+ 


SVM 
Exp. 


-8.380 
-8.481 


1.698 
1.57 


40 


4 (^He) 


(0,0)0+ 


SVM 
Exp. 


-29.937 
-28.295 


1.41 
1.47 


60 


5 


(1,1/2)3/2- 


SVM 
Exp. 


unbound 
unbound 






6 (^He) 
(6Li) 


(n olO+ 
(0,1)1+ 


SVM 
Exp. 
SVM 
Exp. 


-30.07 
-29.271 
-34.59 
-31.995 


2.44 

2.22 
2.43 


600 
600 


TABLE VII. 
H of Eq. (63). 
fm2). 


Energies and rms 
The a-particle is 


radii of N-a systems interactini 
considered to be structureless bi 


I via the Ali-Bodmer potential 
oson. {h^/Ma, = 41.467/4 MeV 


N 


Method 


E (MeV) 




(r2)i/2 (fm) 


K 


2 




unbound 








3 


ATMS g 
SVM 


-5.18 
-5.18 




2.43 
2.43 


60 


4 


ATMS @ 
SVM 


-11.1 
-11.07 




2.65 
2.65 


150 


5 


SVM 


-16.22 




2.99 


400 


6 


SVM 


-20.13 




3.32 


600 
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TABLE VIII. Energies and rms radii of electron-positron systems treated as fermions (f) and 
as bosons (b). Atomic units are used. 



System 




Method 




E 


(r2)i/2 




(e+,e-) 


b,f 


SVM 
Exact 




-0.25 
-0.25 


1.732 
1.732 


10 


(2e+,e- 


) b,f 


SVM 

Variational 
Faddeev § 


41 


-0.262004 
1 -0.2620050702325 
-0.26202 


4.592 
4.594 


150 
700 




h f 


SVM 

kJ V IVl 

Variational 


m 


—0 ^1 ^QSQ 

,|§ -0.515980 


3 608 
3.600 


3on 
300 


(3e+,2e- 
(3e+,2e- 


-)f 
-)b 


SVM 
SVM 




unbound 
-0.5493 


3.53 


1000 
200 


(3e+,3e' 
(3e+,3e^ 




SVM 
SVM 

Variational 


47 


unbound 
-0.820 

1 -0.789 


3.42 


1000 
300 
5 



TABLE IX. Energies and rms radii of "self-gravitating" m-particle-n-antiparticle systems 
(m-f,n— ); f: fermions; b: bosons. VLB and VUB stand for the variational lower and upper 



bounds given 


in Ref. |49|. The 


units of the 


energy 


and length are G'^rrv'h' 


-2 and G- 




respectively. 














System 


Method 




E 






K 


(+,-) b,f 


SVM 




-0.25 




1.732 


10 




Exact 




-0.25 




1.732 




(2+,-) b,f 


SVM 




-1.072 




1.304 


15 




Variational 


49 1 


-1.067 








(2+,2-) b,f 


SVM 




-2.791 




1.027 


100 




VUB (VLB) 


-1.951 


(-3.00) 






(3+,2-) f 


SVM 




-3.758 




1.554 


200 


(3+,2-) b 


SVM 




-5.732 




0.844 


200 




VUB (VLB) 


-4.336 


(-6.25) 






(3+,3-) f 


SVM 




-6.409 




1.621 


300 


(3+,3-) b 


SVM 




-10.215 




0.718 


300 




VUB (VLB) 


-8.130 


(-11.25) 
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Tabic X. Components of the two-particle interaction matrix elements. The operators 
and Cf^ are defined in Eq. (70) of Appendix B. The symbol p specifies the component of 
the two-body interaction. The following additional notations are introduced: ri2 = ri — r2, 
Pi2 = ^(pi - P2), L = ri2 X P12, S = I {ai+ 0-2), 5*12 = 3(<7i • 1-12) (<T2 ■ ri2)- cry (T2, 
X = s'l - s'2, and (x)^ = y^a;Yi^(x). 



p 


definition 


BP 




rip 




c 


1 





\/47r 








t 


S12 













b 


L S 








r -1(1) 
— -i/^i/ XX ( cri+ (T2) 

L J m 





Q 


L2 











<? 




1/2 


V47r!^(-3 + i/x2) 







bb 


(L-S)2 


-i!^2|x2(2+ o-l- 0-2)- 
(x-<Tl) (x- IT2)| 





l/|(x)r„(2+CTl- 0-2)+ 

-ij X X (0-1+ 0-2) 

L J m 

-|(Tlm(x-(T2) + (T2m(x-(Tl))| 


yifL.2|[xXx]S^ + 

r i(2)~i 

2 [xXO-l](l) X [xXO-2](l) \ 

L J J 
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Figure Captions 



Fig. 1 : Different sets of relative coordinates for a system of six identical particles. 



Fig. 2 : Convergence of the ^Li energy on different random paths. The Volkov potential ||3p 
is used. 

Fig. 3 : Convergence of the a-particle energy for the potentials listed in Table I. 



39 



-29.0 
-29.5 
-30.0 
-30.5 
-31.0 
-31.5 
-32.0 



Minnesota 

\ ~ - : 

V ^ ^ - - ATS3_ 

"Vblkov 



MTV 



— \ — \ — I — \ — \ — \ — I — \ — \ — \ — I — \ — \ — \ — I — \ — \ — \ — I — \ — \ — \ — I — \ — \ — \ — I — \ — \ — \ — I 

20 40 60 80 100 120 140 160 



dimension of the basis 




dimension of the basis 



